"crv Development of a Geant4 Solid for Stereo Mini-jet Cells 

^ in a Cylindrical Drift Chamber 

Kotoyo Hoshina a , Keisuke Fujii al , and Osamu Nitoh b , 

> 

a High Energy Accelerator Research Organization(KEK), Tsukuba, 305-0801, Japan 
q ■ fe Tokyo University of Agriculture and Technology, Tokyo 184-8588, Japan 

m ' 
o 

CO 

o 



X 



Oh: 



X 

5-H 



Abstract 



Stereo mini-jet cells will be indispensable components of a future e + e~ linear 
collider central tracker such as JLC-CDC. There is, however, no official Geant4 solid 
available at present to describe such geometrical objects, which had been a major 
obstacle for us to develop a full Geant4-based simulator with stereo cells built in. We 
have thus extended Geant4 to include a new solid (TwistedTubs), which consists of 
three kinds of surfaces: two end planes, inner and outer hyperboloidal surfaces, and 
two so-called twisted surfaces that make slant and twisted ^-boundaries. Design 
philosophy and its realization in the Geant4 framework are described together with 
algorithmic details. We have implemented stereo cells with the new solid, and 
tested them using geantinos and Pythia events (e + e~ — ► ZH at \fs = 350 GeV). 
The performance was found reasonable: the stereo cells consumed only 25% more 
CPU time than ordinary axial cells. 
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1 Introduction 



Experiments at a future linear e + e~ collider such as JLCpP will open up a novel possibility 
to reconstruct all the final states in terms of fundamental particles (leptons, quarks, and 
gauge bosons). This involves identification of heavy unstable particles such as W, Z, and 
t through jet invariant-mass measurements. High resolution energy flow measurements 
will thus be crucial, necessitating high resolution tracking and calorimetry as well as 
good track-cluster matching to avoid double counting. A large cylindrical drift chamber 
with small jet cells (JLC-CDC[2J is our choice for a candidate central tracking device 
to fulfill these requirements. Good track-cluster matching requires, however, small track 
extrapolation errors, which in turn demand high r-0 and z resolutions. The latter forces 
us to introduce cells consisting of stereo wires, since the z resolution in charge division or 
time difference readout is typically 1 % of the wire length or worse, which is a few cm in 
the linear collider use. 

We have already studied and published hardware aspects of common problems in 
designing stereo- wire geometry for a long cylindrical drift chamber with small jet cells . 
In order to finalize the chamber design so as to achieve the best attainable energy flow 
resolution, however, we need to carefully optimize the layout of axial and stereo cells 
through detailed Monte Carlo simulations. Considering the recent advance of object- 
oriented technology in high energy physics software development, we have thus started 
the development of a full Monte Carlo simulator called JUPITERpU based on Geant4[H]. 

The axial layers were easy to implement, using a standard Geant4 solid called 
G4Tubs, which is a segment of a cylinder. On the other hand, there is currently no 
official Geant4 solid available to describe stereo cells, which had been a major obstacle for 
us to install stereo cells into the full simulator. We have thus extended Geant4 to include 
a new solid (TwistedTubs), which comprises three kinds of bounding surfaces: two end 
planes, inner and outer hyperboloidal surfaces, and two so-called twisted surfaces that 
make slant and twisted <fi boundaries. Although TwistedTubs was developed under the 
JUPITER environment, it is actually a general purpose Geant4 extension. This paper 
describes the design philosophy of TwistedTubs, its realization in the Geant4 framework, 
and algorithmic details, together with results of its performance test. 

The paper is organized as follows. We begin with a brief account of geometrical 
parameters that determine the configuration of a stereo cell, and then review the basic 
procedure to add a new solid to Geant4. The subsequent two sections are devoted to de- 
scriptions of design philosophy, implementations, and algorithmic details of TwistedTubs, 
which is followed by presentation of test results and discussions. Finally section |S] sum- 
marizes our achievement and concludes this paper. 

2 Geometry of Stereo Mini-jet Cell 

In this section we introduce, following the convention of 0], the stereo-geometrical pa- 
rameters that will be needed in subsequent sections. 

Consider a cylindrical tube consisting of two cylindrical layers of axial field-shaping 
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Figure 1: An exaggerated illustration of 
a stereo cell as formed by twisting an axial 
cell. 



Figure 2: 3-dimensional view of a single 
stereo wire. 



wires strung across two disc-shaped end plates at some radii, pi n and p out . An axial 
mini-jet cell is a segment cut out from this cylindrical tube by boundaries formed by 
cathode wires. As illustrated in Fig. ^ twisting one end plate by a twist angle A0 with 
the other end fixed turns the axial mini-jet cell into a stereo mini-jet cell. Inspection of 
the figure tells us that the stereo cell consists of three pairs of three kinds of surfaces: a 
pair of hyperboloidal surfaces setting inner and outer radial boundaries, another pair of 
so-called twisted surfaces making left and right azimuthal boundaries, and yet another 
pair of fan-shaped flat surfaces closing the positive and negative z ends of the cell. 

By construction, any of the four side walls of the stereo cell can be regarded as a 
locus of a stereo wire sweeping through the surface. The geometry of a given side surface 
can thus be completely determined by the equation for a single representative wire chosen 
from the stereo wires forming that surface. As depicted in Fig. El the representative stereo 
wire is uniquely specified by the radius at the ends p(z = ±L/2) or that at the center 
p c = p{z = 0), the projected wire length L to the chamber axis, and the twist angle A<p. 
Notice that A0 is signed and measured from A to A\. The stereo angle a, which is 
defined as an angle between AAq and AA\, is also signed, having the same sign as A<p. 

The stereo angle a can now be written in terms of p(z = ±L/2), A<p, and the 
projected wire length (L): 



It is obvious from Fig. El that both the azimuthal angle and the radial position of the 
stereo wire become dependent: 
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p(z) = J {p{z = 0)) 2 + (ztana) 2 

= y/{p{z = ±L/2) f + (z 2 - (L/2) 2 ) tan 2 a, (2.3) 

where z is measured from the middle of the chamber along the chamber axis as Fig. |3] 
indicates. 

Now we consider a hyperboloidal surface obtained as a locus of the straight line 
given by Eqs. 12.21 and 12.31 by sweeping (p(z = 0) from <pi e f t to 4> r i g ht with p{z = 0) fixed. 
Let the outward normal to the hyperboloidal surface at a surface point x = (x, y, z) be n. 
Apparently the outward normal n is in the p-z plane containing x and is perpendicular to 
the tangential vector thereat. Eq. 12.31 tells us that the cross section of the hyperboloidal 
surface by the p-z plane becomes a hyperbola given by 

J (p(z = 0)) 2 + (ztana) 2 
z. 

The tangential vector we need is thus obtained by differentiating this equation with respect 
to z. Normalizing it to |x|, we obtain 

the tangential vector = (ztan 2 a, p), (2.4) 

which leads to the outward normal in the p-z plane: n p = (p, —zt&n 2 a). Recalling that 
p in the p-z plane corresponds to (x, y) in the x-y plane, we now obtain 

n = (x, y, —z tan 2 a). (2.5) 

The inward normal is anti-parallel with this and can be obtained by simply changing the 
signs of all the components. 

On the contrary to the hyperboloidal surface, a twisted surface is a locus of the 
straight line given by Eqs. 12.21 and 12. 3[ when p c = p(z = 0) is swept from p Cjin = 

\J p 2 in - ((L/2)tana) 2 to p C)OUt = \J p 2 out - ((L/2) tana) 2 with <f>(z = 0) fixed, where p in 
and p ou t are the inner and outer radii at the end planes as defined before. By construction, 
any point on such a twisted surface can be specified by two parameters (p c and z) in a 
local coordinate system for which x-axis is chosen such that (f>(z = 0) = 0, as depicted in 
Fig.H 

= Pc 

= p c tan0' = p c Kz (2.6) 
= z, 

where tan0' is given, through Eq. 12. 2\ by 

2z ( A0\ . . 

tan0' = ^■ tan ( _ ^r) ( 2 - 7 ) 
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while k is defined to be 

K = Z' tan (^)' (2 - 8) 

Notice that the A0 must be smaller than ir, implying that p c > in this coordinate 
system. Eq. 12.61 shows that, in addition to the lines corresponding to cathode wires, 
the twisted surface contains another kind of straight lines that appear as the sections of 
different ^-slices of the twisted surface. Each of these lines corresponds to a boundary 
of a z-slice of a stereo cell and passes through the z-axis. At the point x on the surface, a 
pair of these two kinds of straight lines passing through it divides the tangential directions 
at that point into four regions: the surface turns from convex to concave or vice versa as 
crossing the boundaries. In other words, the twisted surface has a saddle shape 2 . 

The two kinds of straight lines contained in the twisted surface are useful when we 
construct a normal to the twisted surface as we shall see below. In general a normal to 
the surface at a given surface point x = (x, y, z) can be formed as a vector product of two 
linearly independent tangential vectors at that point. In our case such tangential vectors 

2 It is well known that the hyperboloidal surface has also a saddle shape. In this case, the two straight 
lines over which the curvature changes its sign are in the directions of stereo angles of ±a. 
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can be readily obtained by partially differentiating Eq. 12.61 with respect to x = p c and z: 





(2.9) 



The first of the above two tangential vectors, e 2 , can virtually be identified as one of the 
cathode wires forming the twisted surface. On the other hand, the second one, e x can 
be viewed as a boundary of a z-slice of a stereo cell. Using Eq. 12.91 we can now easily 
form the surface normal at the surface point x by taking the vector product of these two 
tangential vectors. 



Armed with the equations given in this section, we shall, in the next two sections, 
prepare classes to represent the three kinds of surfaces: J4HyperboloidalSurf aces, 
J4TwistedSurf aces, and J4FlatSurf aces, and assemble them to form a Geant4 solid 
class (J4TwistedTubs) to describe our stereo mini-jet cells of JLC-CDC. 



3 How to Add a New Solid to Geant4 

According to Geant4 User's Guide for Toolkit Developers |5 , every Geant4 solid has to 
inherit from a base class called G4VSolid. This base class has the following pure virtual 
functions: 

G4VSolid 

G4double DistanceToIn (const G4ThreeVector &p) 

to calculate the minimal (or shorter) distance to the solid from an outside point (p), 

G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) 

to calculate the exact distance from the outside point p in the direction of a velocity 
vector v to the solid, 

G4double DistanceToOut (const G4ThreeVector &p) 

to calculate the minimal (or shorter) distance to the solid from an inside point p, 

G4double DistanceToOut (const G4ThreeVector &p, const G4ThreeVector &v, 
const G4bool calcNorm=FALSE , G4bool *validNorm=0, G4ThreeVector *n) 

to calculate the exact distance from the inside point p to the solid along a velocity 
vector v, 

G4ThreeVector Surf aceNormal (const G4ThreeVector &p) 

which returns the outward unit normal at a surface point p (or, if p is not on the 
surface, at the surface point that is the closest from p), 
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EInside Inside (const G4ThreeVector &p) 

to judge whether a space point p is inside, or outside, or on the surface of the solid 
and returns klnside, or kOutside, or kSurface, accordingly, 

G4bool CalculateExtent (const EAxis pAxis, const G4VoxelLimits fepVoxelLimit , 

const G4Af f ineTransf orm fepTransform, G4double fepMin, G4double fepMax) const 
to calculate minimum and maximum extents (pMin and pMax) of the solid in the di- 
rection of a given coordinate axis pAxis within the limits specified by pVoxelLimit 
under an Affine transformation pTransf orm, and 

G4GeometryType GetEntityType () const 

to identify geometry type of the solid (necessary for persistency and STEP interface, 
but otherwise unused). 

The first five involve surfaces that bound the solid and can be calculated on a surface-by- 
surface basis. On the other hand the remaining three require information of the solid as a 
whole. Any new user defined solid thus has to be equipped with these functions. Geant4 
also provides BREPS 3 classes for the purpose of constructing a solid from boundary 
surfaces. BREPS are, however, primarily designed to facilitate implementation of a solid 
with a complicated shape through interface like STEP to a CAD system. Since the stereo 
cell geometry can be handled analytically as explained in the subsequent sections, we 
decided to develop a new dedicated solid class for it. 

4 Design of J4TwistedTubs 

Before coding our new solid " J4TwistedTubs", we set the following guide line to fulfill 
the required functionality discussed in the last section: 

1. J4TwistedTubs should be implemented as a collection of bounding surfaces. It will 
thus have to have, as its data member, an array of pointers to the instances of the 
corresponding surface classes. 

2. All of these surface classes must inherit from an abstract class named " J4VSurf ace" , 
which carries generic information on each surface and sets basic interfaces to all the 
surface classes that inherit from it. 

3. Each surface is orientable and should know which side is outside, but the calculation 
of the distance to the surface from a point p should not depend on whether the point 
is inside or outside of the solid, since the distance should only depend on the shape 
of the surface in question. 

4. It is then the role of "DistanceToIn" or "DistanceToOut" methods of its parent 
J4VSurf ace class to judge the distance depending on the context such as the angle 
between the surface normal and the velocity vector v. This way, we can avoid 
repetition of the same code. 

3 Boundary REPrcscntcd Solid. 



6 



J4TwistedTubs 



DistanceToIn(p, v) 

DistanceToIn(p) 

DistanceToOut(p, v) 

DistanceToOut(p) 

SurfaceNormal(p) 

Inside(p) 



J4VSurface* fSurface[6] 



J4FlatSurface 



DistanceToSurface(p, v) 
DistanceToSurface(p) 
GetNormal(p) 
GetAreaCode(p) 



G4Three Vector fNormal 



J4VSurface 



DistanceToIn(p, v) 

DistanceToOut(p, v) 

DistanceTo(p) 

SurfaceNormal(p) 

virtual DistaceToSurface(p, v) =0 

virtual DistaceToSiuface(p) =0 



J4VSurface *£Neighbours[4] 
G4int fHandedness 
G4Three Vector fTranslate 
G4RotationMatrix fRot 



n n 



J4TwistedSurface 



DistanceToSurface(p, v) 
DistanceToSurface(p) 
GetNormal(p) 
GetAreaCode(p) 



G4double fKappa 



J4HyperbolicSurface 



DistanceToSurface(p, v) 
DistanceToSurface(p) 
GetNormal(p) 
GetAreaCodefp) 



G4double fKappa 
G4double fRO 



Figure 5: A class diagram in UML of 
J4TwistedTubs 



f Orientation = -1 



fOrientation = 1 




fOrientation = 1 



fOrientation = 1 



fOrientation = -1 



fOrientation = -1 



Figure 6: Values of 

fOrientation data member 
of J4TwistedTubs 



5. All J4TwistedTubs has to do will then be to simply invoke "DistanceToIn" or 
"DistanceToOut" for each of the bounding surfaces and choose the best 4 . The sur- 
face normal can then be calculated by calling the corresponding GetNormal method 
for the selected surface. 

As described in section EJ J4TwistedTubs consists of three pairs of three kinds of sur- 
faces. We named these three kinds of bounding surfaces J4HyperboloidalSurf ace, 
J4TwistedSurf ace, and J4FlatSurf ace, all of which are descendants of J4VSurface. 
Interrelation of these classes is illustrated in Fig. as a class diagram. On the other 
hand, Fig. El shows the outward direction of each surface as stored in a data member 
called fOrientation. 

As sketched above, the abstract base class J4VSurf ace has the following methods: 

J4VSurface 

G4int DistanceToSurf ace (const G4ThreeVector &p, G4ThreeVector *xx, 
G4double *distance, G4int *areacode) = 

which calculates the distance (*distance) from a point p to the surface as well as 
the point (*xx) of the closest approach to the surface and an *areacode given by 
GetAreaCode explained below, and then returns the number of the points of closest 
approach that is always 1. 

4 The best is usually the one that returned the shortest distance. One exception is the case in which 
a particle is coming into the solid from a corner or edge of the solid. 
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G4int DistanceToSurf ace (const G4ThreeVector &p, const G4ThreeVector &v, 
G4ThreeVector xx[] , G4double distance [], G4int areacode[], 
G4bool isvalid[] , EValidate validate) = 

which calculates the distance(s) (distance []) along the velocity vector v to the 
surface from the point p, the intersection(s) (xx [] ) of the particle track with the 
surface, and the areacode(s) (areacodeG) of intersection(s), and then returns the 
number of intersections. If so flagged by validate, it does boundary check for each 
of the candidate intersection(s) and passes the test result through isvalidG. The 
return value then becomes the number of valid intersection(s). 

G4ThreeVector GetNormal (const G4ThreeVector fax, G4bool isglobal) = 

which calculates and returns the normal at a surface point xx. The surface point 
xx and the normal to be returned are assumed to be given in the coordinate system 
of its mother solid or in the local frame of the surface in question, depending on 
whether isglobal is TRUE or not, respectively. 

G4Int GetAreaCode (const G4ThreeVector fax, G4bool withTol) = 

which decides whether the intersection xx given by DistanceToSurf ace (p,v) is 
on a boundary, or at a corner, or inside or outside of the surface with or without 
tolerance as flagged by withTol, and then returns a corresponding area code. 

G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v, 
G4ThreeVector fax) 

which invokes DistanceToSurf ace (p,v) of its descendant concrete class and judges 
the validity of each of the resultant intersection(s) and the corresponding distance(s), 
examining the sign of the distance, the angle of the particle velocity v to the surface 
normal at the intersection point on the surface. If the intersection is on a boundary 
or at a corner of the surface, it also checks the angle to the surface normal (s) of 
the adjacent surface(s). It then returns the distance if valid, or infinity otherwise, 
together with the valid intersection xx, if any. 

G4double DistanceToOut (const G4ThreeVector &p, const G4ThreeVector &v, 
G4ThreeVector fax) 

which invokes DistanceToSurf ace (p) of its descendant concrete class and judges 
the validity of each of the resultant intersections(s) and the corresponding dis- 
tance^), taking into account the sign of the distance and the angle between the 
surface normal and the particle velocity v, and then returns the distance to the 
caller if valid, or infinity otherwise, with the valid intersection passed through xx, 
if any. 

G4double DistanceTo (const G4ThreeVector &p) 

which just passes the return value of DistanceToSurf ace (p) without any additional 
judgment. This functions is used in DistanceToIn(p) or DistanceToOut (p) of 
J4TwistedTubs. 
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Notice that the first four of these member functions are pure virtual and should 
be implemented in its derived concrete classes. Since we have already explained how to 
construct a surface normal in section |21 and since description of GetAreaCode would only 
involve technical details, we will concentrate, in what follows, on how we implemented 
the base class J4VSurf ace and the first two that are the DistanceToSurf ace functions 
with and without the particle velocity, for each of the three derived surface classes: 
J4HyperboloidalSurf ace, J4TwistedSurf ace, and J4FlatSurf ace. 

5 Algorithm and Implementation 

5. a J4VSurface 

Being the base class for all the surface classes that make J4TwistedTubs, the J4VSurf ace 
class plays the following four major roles: (1) it carries data members to store generic 
information on a surface, (2) it standardizes interface for generic methods to be imple- 
mented by its descendant, (3) it provides tools commonly used by any derived surface 
class, and (4) it steers the distance calculation in DistanceToIn and DistanceToOut, by 
filtering the return values from DistanceToSurf ace of its derived class. We now describe 
these points below in some more detail. 

5.a.l Generic Data Members 

J4VSurf ace has the following protected data members to store generic information on its 
derived surface class. 

• Two axes of the surface 

• Maximum and minimum limits of the surface along the two axes 

• The rotation matrix and the translation vector to transform the surface from its 
mother solid (J4TwistedTubs) coordinate system to the local frame attached to the 
surface 

• The last value of the distance to the surface from a point p 

• The last value of the surface normal 

On the other hand, the next five data members are private, since they are basically used 
internally by J4VSurf ace: 

• Four pointers to neighboring surfaces 

• Four positions of the corner points 

• An instance of a local class that stores information on the four boundary line seg- 
ments 
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• The name of the surface 



• A pointer to its mother solid (J4TwistedTubs) 

Notice that J4VSurf ace currently assumes the existence of four corners and consequently 
four boundary line segments for its derived surface class. It is, however, relatively easy to 
extend it to support variable number of neighboring surfaces. 

5.a.2 Common Utility Functions 

G4double DistanceToLine (const G4ThreeVector &p, const G4ThreeVector faO, 
const G4ThreeVector fed, G4ThreeVector fax) 

to calculate the distance from a point p to a line specified by a reference point xO on 
it and the line direction vector d, together with the point xx of the closest approach 
on the line. It then returns the distance. 

G4double DistanceToPlane (const G4ThreeVector &p, const G4ThreeVector &x0, 
const G4ThreeVector &n0, G4ThreeVector fax) 

to calculate the distance from p to a plane specified by a normal nO and a reference 
point xO on it, together with the point xx of the closest approach to the surface. It 
then returns the distance. 

G4double DistanceToPlane (const G4ThreeVector &p, const G4ThreeVector faO, 
const G4ThreeVector &tl, const G4ThreeVector &t2, G4ThreeVector fax, 
G4ThreeVector &n) 

to calculate the distance from p to a plane and the point xx of the closest approach 
to the plane, where the plane is specified by a reference point xO and two linearly 
independent vectors tl and t2 on the plane. It then returns the distance. 

G4double DistanceToBoundary (G4int areacode, G4ThreeVector fax, 
const G4ThreeVector &p) 

to calculate the distance from p to a boundary specified by areacode, and the point 
xx of the closest approach to the boundary. Its return value is the distance. 

G4int AmlOnLef tSide (const G4ThreeVector feme, const G4ThreeVector &ref , 
G4bool withTol) 

to judge relation between two vectors me and ref by comparing their angles. If 
me lies clearly on the left side of ref, namely me has a angle smaller than that 
of ref, it returns 1. If withTol is TRUE and if me lies within ±0.5xkAngTolerance 
(= 0.5 x 10~ 9 ) radians of ref, it returns 0. Otherwise, it returns -1. 

5.a.3 DistanceToIn and DistanceToOut 

DistanceToIn(p, v) 

On invocation, DistanceToIn (p,v) 5 proceeds according to the following algorithm: 
5 In what follows, we often leave out arguments which are irrelevant in the explanation of algorithmic 
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1. Initialize a temporary data buffer to store the current best distance and intersection 
to klnf inity. 

2. Then call DistanceToSurf ace(p,v) and get a set of distance(s) (D) and intersec- 
tion^) (X): the number of intersections is at most 2 in our case. 

3. If a candidate D and its corresponding X satisfy the following conditions, discard 
this candidate and move on to the next, if any: 

• D is negative, or 

• the scalar product between the outward surface normal at X and the particle 
velocity v is positive, or 

• X is outside of the surface. 

4. If X is inside of the surface and if D is smaller than the last value stored in the 
temporary data buffer defined in step 1, update the stored D and X and move on 
to the next intersection, if any. 

5. If X is on the boundary of the surface within ±0.5xkCarTolerance(= 0.5 x 10~ 9 ), 
we need to judge whether the particle is really entering the solid or just grazing the 
solid and flying away. In this case, invoke the DistanceToSurf ace (p,v) method of 
the (two) neighboring surface (s) (two if the intersection is at a corner, one otherwise) 
and get a set of distance(s) (D f ) and intersection(s) (X'), and if X' satisfies the 
following conditions, exit DistanceToIn immediately by returning klnf inity as 
the distance to the surface originally in question. 

• If X' is inside of a neighboring surface 6 , or 

• if the X' is on the same boundary (or at the same corner) as X and if the 
scalar product of the outward surface normal at X' and the particle velocity v 
is positive 7 . 

6. If X survives all the tests for all the neighboring surfaces, check whether the current 
D is smaller than that stored in the temporary data buffer, and if so, update the 
stored D and X. 

7. Return the stored distance in the temporary data buffer. 
The flow chart of DistanceToIn is shown in Fig. 

DistanceToOut (p , v) 

On the other hand, the algorithm for DistanceToOut (p,v) is less complicated as 
described below: 

aspects of each function, as long as there is no possibility of confusion. 

6 In this case this neighboring surface should be the one through which the particle enters the solid. 

7 A particle cannot come into a solid from the edge (or corner) if the scalar product is positive for any 
of the surfaces forming the edge (or corner). 
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i = 0; best = klnfinity; 
nx = DistanceToSurface 

D[nx],X[nx], etc. 

(nx: number of X) 



best = D[i] 




Loop B 



return best 



B 



j = 0; 

nsurf = GetNeighborSurfacesO 
(nsurf: number of neighbour 
surfaces at present boundary 
1: boundary, 2: corner) 



k = 0; 

nx' = DistanceToSurface 
D'[nx], X'[nx], etc. 
(nx' : number of X') 



k++ 




return klnfinity 



Figure 7: Flow chart of Distance Toln 
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1. Initialize a temporary data buffer to store the current best distance and intersection 
to klnf inity. 

2. Call DistanceToSurf ace(p,v) and get a set of distance(s) (D) and intersection(s) 
(X): the number of intersections is at most 2 in our case. 

3. If D and X meet the following conditions, reject this candidate and move on to the 
next: 

• D is negative, or 

• the scalar product of the outward surface normal at X and the particle velocity 
v is negative, or 

• X is outside of the surface. 

4. If D is smaller than the stored value in the temporary buffer, update the stored D 
and X. 

5. Return the best distance stored in the temporary buffer. 
DistanceTo(p) 



DistanceTo(p), which is called by the DistanceToIn(p) and DistaceToOut (p) meth- 
ods of J4TwistedTubs, just returns DistanceToSurf ace (p) which is implemented in the 
derived surface class. 

5.b J4HyperboloidalSurface 

Geant4 provides a solid class called G4Hype to represent a hyperboloidal volume with 
its symmetry axis parallel to the z-axis. Since J4HyperboloidalSurf ace follows the 
basic algorithm to calculate the inward or outward normal and the minimal distance to a 
hyperboloidal surface, we only sketch it here. 

5.b.l DistanceToSurface(p, v) 

The distance to a hyperboloidal surface along a velocity v = (v x , v y , v z ) from a point 
P = {PxiPyiPz) can be obtained as follows. Define the intersection of a particle track with 
the hyperboloidal surface to be X = p + £v = (p x + tv xi p y + tv y ,p z + tv z ), where t, the 
time, equals the distance we want, provided that v is normalized to unity. From Eq. 12.21 
we then have a quadratic equation for t: 

ip x + tv x ) 2 + (p y + tVy) 2 = r 2 z=0 + (p z + v z ) 2 tan 2 a. (5.10) 
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p=z tana 



dp/dz = ztan 2 a/p 




From Inside 



Pi 



O 



X7p 



Figure 8: Distance from outside 



Figure 9: distance from inside 



5.b.2 DistanceToSurface(p) 

DistanceToSurf ace(p) is to calculate the minimum distance to the surface in question. 
Since this function is primarily used to set the safety radius not to cross any volume 
boundaries, it is allowed to underestimate the minimum distance. In the calculation of 
this minimum distance we can thus approximate the surface curve by a properly chosen 
straight line. Taking advantage of z-symmetry, we only have to consider the region of 
z > 0. 

Distance from Outside 

In this case, the proper line can be chosen in the following way: 

1. Find the surface point (Xi) which has the same z as p. 

2. Draw a line from p which is perpendicular to the asymptotic line p = z tana and 
let their intersection be Q. 

3. Find the surface point (X 2 ) that has the same z as Q. 

4. Connect Xi and X 2 

Then, by construction, the distance from p to the line X X X 2 never exceeds the mini- 
mum distance to the hyperboloidal surface. The z component of X 2 can be obtained 
geometrically from Fig. |H1 



Distance from Inside 

In this case, a proper line will be the tangential line at the surface point Xi which has 
the same z component as p (see Fig. |UJ). The tangential vector is given by Eq. 12.41 



z 2 = OQ cos a 
OQ = OR cos a 
OR = zi+pp tana. 



(5.11) 
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5.c J4FlatSurface 

5.C.1 DistanceToSurface(p, v) 

The distance to a plane along a velocity vector v from a point p is readily obtained from 



where n is the normal to the plane and xo is the reference point on the plane. 
5.C.2 DistanceToSurface(p) 

Since the two planar boundary surfaces of J4TwistedTubs are perpendicular to the z- 
axis, which means its normal is in the direction of the z-axis. All DistanceToSurf ace(p) 
has to do is simply to return the difference between z-components of p and the plane in 
question. 

5.d J4TwistedSurface 
5.d.l DistanceToSurface(p, v) 

The distance to a twisted surface from a point p = (p x ,P y ,Pz) along a velocity vector 
v = (v x ,v y ,v z ) can be obtained in a similar manner to that for J4HyperboloidSurf ace. 
Substituting X = (X x , X y , X z ) = p + tv in Eq. 12. 6^ we have the following quadratic 
equation: 



for the distance t to the surface. 

Notice that, although this quadratic equation may contain a solution for which 
X x = p c < 0, such a situation should not take place in practice, since the twist angle 
cannot exceed A0 = tc without breaking wires. This requires X x = p c in Eq. 12.61 be 
positive. Notice also that, for a given stereo angle, the twist angle A0 increases with the 
wire length. In solving Eq. 15.131 we thus need to set limits on z as determined by the 
endcap locations. These validity checks can be made for the obtained crossing point X 
after solving the above equation. 

5.d.2 DistanceToSurface(p) 

The exact equation for the minimal distance to the twisted surface from a point p is 
biquadratic. As in the case of J4HyperboloidalSurf ace, what we really need is, however, 
not the exact minimal distance. Its approximation will do, provided that it never exceeds 
the true minimal distance. We thus developed an approximation method that uses the 
distance to an appropriately chosen plane instead of solving the exact equation. Now 
the question is how we find such a plane to approximate the twisted surface. Fig. [TU] 
illustrates the procedure to find it. 

Since we should never overestimate the distance, the plane has to lie on the same 
side as the point p with respect to the twisted surface, at least in the vicinity of the 






(5.13) 
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Figure 10: Procedure to choose an appropriate plane to approximate the minimum 
distance to the twisted surface 
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point of the closest approach to the plane. In order to fulfill this condition, essential is 
the following observation about the geometric features of the twisted surface. As stressed 
in section |21 the twisted surface contains two kinds of straight lines given by e 2 and e x 
in Eq. 12.91 one that can be regarded as one of the field wires making up the twisted 
surface and the other that can be taken as a boundary of a z-slice of a stereo cell. Lines 
AD and BC of Fig. El are the first type (e 2 -type) , while AB and DC are the second 
type (e x -type). We can thus cut out, from the twisted surface, a segment ABCD that is 
bounded by four straight lines. Recall that this segment has a saddle shape: it is concave 
along the diagonal AC in Fig. while it is convex along the other diagonal BD 8 . 

Our goal is to choose a plane that is placed in front of both the current point p and 
the twisted surface and to calculate the distance to that plane. Once we find a diagonal 
which does not cross the twisted surface, we can immediately span a plane with the 
diagonal and an adjacent side (AB or DC) which is entirely lies on the surface 9 . Hence 
our first task is to select an appropriate diagonal along which the surface is concave. The 
procedure in DistanceToSurf ace(p) to select such a diagonal can be itemized as follows: 

1. Let A and C be the points of the closest approach from the point p to the inner 
and outer boundary lines (m and n in Fig. ITU1-I). If the z components of A and C 
do not sandwich the z component of p, move A or C along the boundary lines until 
its z component matches that of p. 

2. Choose D and B on lines m and n, respectively, in such a way that lines AB and DC 
are perpendicular to the z-axis and thus entirely contained in the twisted surface. 

3. First check if diagonal AC satisfies our requirement by examining the sign of product 
of the following three test variables: 

(a) the signed twist angle A<p, 

(b) an orientation index, which is 1 if p is on the negative <fi side of the twisted 
surface, if p is on the surface within ±0.5xkAngTolerance (= 0.5 x 10 -9 ) 
radians, and -1 otherwise, and 

(c) the z coordinate of A minus that of C. 

If the product is positive, accept the diagonal. If it is 0, which means either (b) 
is or (c) is 0. If (b) is 0, p is on the surface (or is very close to the surface). In 
this case, immediately exit DistanceToSurf ace (p) by returning as the minimal 
distance to the twisted surface. Else if (c) is 0, return the distance to line AC from 
p. Finally if the product turns out to be negative, which means that the diagonal is 
invalid crossing the twisted surface, exchange A and D, and C and B, respectively. 

Diagonal AC can now be used to set up two planes: plane ADC and ABC. The normal 
to these planes are indicated in Fig. HU1-L and should point into the hemisphere that 

8 Of course, which is convex or concave flips, depending on from which side you are looking. 
9 Remember the twisted surface can also be regarded as a locus of line AB in Fig. 0] 
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contains p. The subsequent distance calculation has two branches, depending on the 
relative location of p with respect to these two candidate planes. 

Case A: p is in region i), or ii), or iii) of Fig. fTTH-h 



In these cases, all DistanceToSurf ace(p) needs to do is to return the smallest positive 
distance. 

Case B: p is in region vi) of Fig. [TUl-r 

In this case, we must re-define the planes so that at least one of them lies behind point 
p. Let M and N be the middle points of AB and CD, respectively. Since line MN is 
fully contained in the twisted surface, we can now define new quadrangles ADNM and 
MNCB. We can then define four planes ADN, AMN, MNC, and MBC. Planes ADN 
and MBC, however, coincide with ADC and ABC, respectively, and have hence been 
tested already. Therefore, we just have to test the remaining two planes. Again there are 
two cases to consider, depending on the relative location of point p with respect to the 
planes (see Fig. HUlI'). 

Case B-l: p is in region i) of Fig. [TUl-ir 

In this case, DistanceToSurf ace (p) should just return the smaller of the distance 
to AMN and that to MNC. 

Case B-2: p is in region ii) or iii) of Fig. ITUl-ir 

In these cases, the situation is essentially the same as in iv) of Fig. HUH'. Compare 
the distances to AMN and MNC and split into two the one that gave a negative 
distance. 

Notice that, by construction, p must not fall into region iv) of Fig. ITUl-ir. If it happens, 
therefore, we just abort the program. 

The procedure for B-2 is actually implemented as a recursive call: J4TwistedSurf ace is 
equipped with its own DistanceToPlane, which calls itself recursively until the situation 
changes to B-l and returns a positive distance. In this way DistanceToSurf ace can get 
the desired distance by a single call to DistanceToPlane, when case B takes place. 

6 Test of J4TwistedTubs 

In order to test the new solid class J4TwistedTubs, we prepared a test program (for stereo 
cells) as follows. We first constructed a world volume with a cubic shape of 6m x 6m 
x 6m and put, at its center, a cylindrical tube- type layer with inner and outer radii of 
30cm and 130cm, respectively, and a length of 260cm. In this cylindrical layer, we placed 
a J4TwistedTubs object, which has a twist angle (A0) and a 0-width of tt/3, inner and 
outer radii at the endcaps of 50cm and 100cm, respectively, and a z-length of 200cm. Into 
this mother J4TwistedTubs volume, installed were two daughter J4TwistedTubs objects 
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Figure 11: Schematic view of the test geometry 



which have half the 0-width but otherwise the same geometry. All of these volumes are 
made of air. A schematic view of the test geometry is shown in Fig. ^2 

In order to evaluate the computing speed of J4TwistedTubs, we also prepared an- 
other test program (for axial cells) with the mother and the two daughter J4TwistedTubs 
objects replaced by ordinary cylindrical tubs (G4Tubs). 

6. a Test with Geantinos 

Geant4 provides a hypothetical particle called " geantino" designed primarily for geometry 
tests. A geantino does not interact with any materials and flies along a straight line 
trajectory. It only makes a hit when it crosses a boundary of volumes which are sensitive 10 . 
We first carried out a geantino test of J4TwistedTubs, using the aforementioned geometry 
test program with only the daughter J4TwistedTubs objects made sensitive. Fig. 1121 is 
a 3-dimensional view of hit points made by 10000 geantinos injected from 4ir steradian 
around the surrounding cylinder. As seen from the figure all the hit points are correctly 
on the surface of the daughter J4TwistedTubs volumes. 

6.b Test with Higgs Events 

Unlike geantinos, a real charged particle makes a curved trajectory in a central tracker 
due to the magnetic field applied to it, and interacts with various detector materials, 
producing, for instance, low energy tracks such as 5-rays which might complicate the 
particle tracking through the detector volumes. In order to stress-test our J4TwistedTubs 

10 Geant4 generates a hit when (1) a physics process took place on a flying particle or (2) its track 
crossed a boundary of a volume registered as a sensitive detector. 
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Figure 12: 3- dimensional view of hit points of geantinos on the daughter J4TwistedTubs 
volumes 




Figure 13: x-z view Figure 14: x-y view 



under more realistic environment, we also tried Higgs events (e + e~ — > ZH at ^/s = 
350 GeV) generated by Pythia 6 j. In this case, the number of tracks with a transverse 
momentum of 1 GeV or greater is around 50 per event. Because of the lack of a calorimeter 
in this test program, however, relatively low energy tracks curl up and pass through 
the sensitive volumes multiple times. Taking into account this curl-up effect and the 
geometrical acceptance of the J4TwistedTubs volumes being about (7r/3)/(27r) = 1/6, 
the average number of tracks passing through the twisted volumes was estimated to be 
around 10 per event. Processing of 1000 Higgs events, which thus correspond to roughly 
10000 tracks hitting the sensitive volumes, took 11647 seconds on a Power Macintosh 
800MHz (Dual CPU) with 1GB memory for the J4TwistedTubs test program. The same 
test took 9242 seconds on the same platform for the axial cell test program. The 25% 
extra CPU time consumed by J4TwistedTubs is acceptable for our purpose, considering 
the expected CPU time necessary for calorimeter simulation. 

7 Discussion 

So far we treated the twisted surface analytically in J4TwistedSurf ace. We found, how- 
ever, Eq. 15.131 which is the key equation for our analytic treatment, suffers from roundoff 
errors in the situation sketched below. 

In general a roundoff error enters as a result of a subtraction of two numbers that 
share many common digits. In the case of J4TwistedSurf ace, this happens when v x or 
v z or both become small, and consequently the second term of the discriminant of the 
quadratic equation Eq. 15.131 becomes negligible. Geometrically v x = corresponds to an 
extreme case in which the velocity vector of the particle is contained in a plane that is 
spanned by a straight line in the direction of a wire and another straight line which is 
parallel with the z-axis. On the other hand, v 2 = is the case in which the velocity vector 
is in the plane that is perpendicular to the z-axis. In both of these cases, the problem 
becomes purely 2-dimensional and hence the equation becomes linear. 

When the quadratic equation is used in such a case, DistanceToSurf ace(p, v) 
might return a crossing point that falls short of or goes beyond the surface by more than 
kCarTolerance. The next step might then end up with an impossible situation where 
DistanceToIn (DistanceToOut) would be called from inside (outside) of the volume, 
confusing the tracking program in Geant4. 

When the calculated crossing point (X orig in Fig. [T3J) is found to be more than 
kCarTolerance-off the surface, we apply the following Newtonian method in order to 
avoid such inconsistency. We first approximate the twisted surface, at the surface point 
that has the same local x and z coordinates (X sur f), by a tangential plane spanned by 
two straight lines contained in the twisted surface: one in the wire direction and the other 
in the radial direction. We then recalculate the crossing point (X new ) with this tangential 
plane. We iterate this procedure until the recalculated crossing point is found on the 
surface within the tolerance. 
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Figure 15: Procedure to fix the estimated intersection (X orig ) when it falls short of 
the surface due to a roundoff error. The surface point (X sur f) has the same x and z 
coordinates as X orig . We approximate the twisted surface by the tangential plane at 
X sur f, which is spanned by two straight line sections l\ and I2, and calculate a new 
crossing point (X new ) with this plane. In the limit of v z — or v x = 0, for which the 
numerical instability is expected, the new intersection will be on l 2 or li, respectively, and 
therefore the approximation will become exact. 
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8 Summary and Conclusion 



We have developed a new Geant4 solid called J4TwistedTubs in order to handle stereo 
mini-jet cells of a cylindrical drift chamber like JLC-CDC. This new solid consists of three 
kinds of surfaces, each of which is represented respectively by J4HyperboloidalSurf ace, 
J4FlatSurf ace, and J4TwistedSurf ace. These three surface classes correspond to inner 
and outer hyperboloidal surfaces, two end planes, and two so-called twisted surfaces that 
make slant and twisted 0-boundaries, respectively 

There has been no Geant4 object for the twisted surface. In this paper, we have 
thus explained the algorithmic details of our new surface class (J4TwistedSurf ace) for 
it. It should also be emphasized that all of these three surface classes are derived from 
a single base class called J4VSurf ace, which greatly facilitated the distance and surface 
normal calculations. 

We have implemented stereo cells with the new solid, and tested them using geanti- 
nos and Pythia events (e + e~ — > ZH at ^/s = 350 GeV). The stereo cells consumed 25% 
more CPU time than ordinary axial cells did. We found this acceptable, considering the 
expected CPU time necessary for calorimeter simulation. 
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